General Disclaimer 


One or more of the Following Statements may affect this Document 


• This document has been reproduced from the best copy furnished by the 
organizational source. It is being released in the interest of making available as 
much information as possible. 


• This document may contain data, which exceeds the sheet parameters. It was 
furnished in this condition by the organizational source and is the best copy 
available. 


• This document may contain tone-on-tone or color graphs, charts and/or pictures, 
which have been reproduced in black and white. 


• This document is paginated as submitted by the original source. 


• Portions of this document are not fully legible due to the historical nature of some 
of the material. However, it is the best reproduction available from the original 
submission. 


Produced by the NASA Center for Aerospace Information (CASI) 



A major purpose of the Techni- 
cal Information Center is to provide 
the broadest dissemination possi- 
ble of information contained in 
DOE’s Research and Development 
Reports to business, industry, the 
academic community, and federal, 
state and local governments. 

Although a small portion of this 
report is not reproducible, it is 
being made available to expedite 
the availability of information on the 
research discussed herein. 

1 




DOE/JPL/95584 3- -8 3/8 
DE8 3 008C71 


DIST. CATEGORY UC-63 
DOE/JPL-955843/83/8 
DRD No. SE5 
DRD No. 139 

LARGF-AREA SHEET TASK 

ADVANCED DENDRITIC- WEB -GROWTH DEVELOPMENT 

C. S. Duncan, R. G. Seidenst icker , J. P. McHugh, and 
J. Schruben 

Annual Report. 

October 23, 1981 to October 22, 1982 
January 18, 1983 

Contract No. 955843 

The JPL Flat Plate Solar Array Project is sponsored by the U. S. 
Dept, of Energy and forms part of the Solar Photovoltaic 
Conversion Program to initiate a major effort toward the 
development of low-cost solar arrays. This work was performed for 
the Jet Propulsion Laboratory, California Institute of Technology, 
by agreement between NASA and DOE. 


NOTICE 

PORTIONS OP . HIS REPORT ARE ItLECIBlE. * 

It has bein reproduced from the best^ 

available copy to permit the broadest] ' 

possible availability. 


Westinghouse R&D Center 
1310 Beulah Road 
Pittsburgh. Pennsylvania 15235 


DIST. CATEGORY UC-63 
DOE/JPL-955843/83/8 
DRD No. SE3 
DRD No. 139 


LARGE -AREA SHEET TASK 

ADVANCED DENDRITIC WEB GROWTH DEVELOPMENT 


C. S. Duncan, R. G. Seldensticker , J. P. McHugh, and 
J. Schruben 

Annual Report 

October 23, 1981 to October 22, 1982 
January 18, 1933 

Contract No. 953843 

The JPL Flat Plate Solar Array Project Is sponsored by the U. S 
Dept, of Energy and forms part of the Solar Photovoltaic 
Conversion Program to initiate a major effort toward the 
development of low-cost solar arrays. This work was performed for 
the Jet Propulsion Laboratory, California Institute of Technology, 
by agreement between NASA and DOE. 


DISCLAIMER 


Approved : 


Jll&s 


R. Mazelsky" Manager 
Crystal and Device Technology 


'f 


This report was prepared as an account of work spon-orcd by an agency of the United States 
(iovernmcn Neither the United States Government nor any agency thereof, nor any of their 
employees makes any warranty, express or implied, or assumes any legal liability or responsi- 
bility U th-- accuracy, completeness, or usefulness of any information, apparatus, product, or 
process disclc cd, or represents that its use would not infringe privately owned rights. Refer- 
ence herein to any specific commercial product, process, or service hy trade name, trademark, 
. lanufi.cturer, or otherwise docs not necessarily constitute or imply its erjorsement, recom- 
mendation or favoring by the United States Government or any agency thereof The views 
and opinions of authors expressed herein do not necessarily state or reflect those of the 
Untied States Government or any agency thereof 


Approved: 


s — ; — : 

R. T. Begley, Manager 
Solid State, K.&D 


Westingnouse R&D Center 
1310 Beulah Road 
Pittsburgh. Pennsylvania 15235 



CONTENTS 


LIST OF FIGURES iii 

1 . SUMMARY 1 

2. INTRODUCTION 2 

3. PROGRESS IN WEB GROWTH RESEARCH 4 

3.1 Temperature and Stress Modeling 4- 

3.1.1 Introduction ....4 

3.1.2 Model Development 9 

3. 1.2.1 Introduction 9 

3. 1.2. 2 Review of Radiative Heat Transfer.. 10 

3. 1.2. 3 New Geometry 13 

3. 1.2. 4 Verification of the New Model..... 16 

3.1.3 Operation of the Model 18 

3.1.4 Stress Model Development 23 

3.1.5 Synthetic Temperature Profiles 26 

3.1.6 Development of New Configurations.... 32 

3.2 Operation of the New Experimental Web Growth Facility ..41 

3.3 Experimental Web Growth .........44 

3.3.1 Introduction 44 , 

3.3.2 The J419 and J435 Configurations 44 

3.3.3 The J460 Configuration 47 

3.3.4 Long Growth Runs and Oxide Control .49 

4. SUMMARY, CONCLUSIONS, AND FUTURE WORK 50 

4.1 Summary and Conclusions ...50 

4.2 Future Work. 50 

5 . NEW TECHNOLOGY 51 

6. REFERENCES 52 

7 . ACKNOWLEDGEMENTS 53 

8. APPENDICES ...54 


i 


LIST OF FIGURES 


Page 

Figure 1 Application of computer models 8 

Figure 2 Geometrical and thermal parameters of previous model 

for calculating temperature profile in the web 11 

Figure 3 Radiation geometry . . 12 

Figure 4 Shading limits for viewing the funtace elements from 

the web* . . » . * 14 

Figure 5 Viewing regions from a point on the web..... . ... 17 

Figure 6 Model geometry and (aT)" for run 9-11C * 19 

Figure 7 Model geometry and (aT)" for run 9-11CA (beveled 

elements) 20 

Figure 8 Example of numerical output from new web temperature 

computer code 21 

Figure 9 Constant stress temperature profiles 27 

Figure 10 End effect region of constant stress profile 29 

Figure 11 Stresses generated by abruptly changing zero stress 

temperature profiles 31 

Figure 12 Temperature model of J98M3 configuration 34 

Figure 13 Delta x-stress for J98M3 configuration. 35 

Figure 14 J419 configuration 36 

Figure 15 Delta x-stress for J419 configuration 37 

Figure 16 Temperature profile results for J460 configuration 39 

Figure 17 Delta x-stress distribution for J460 configuration 40 

Figure 18 New experimental web growth furnace 42 

iii 



1. SUMMARY 


During this reporting period, substantial progress has been made 
in the reduction of thermally generated stresses in the growing web 
crystal. These stresses, which if too high cause the ribbon to 
degenerate, have been reduced by a factor of three, resulting in the 
demonstrated growth of high-quality web crystals to widths of 5.4 cm. 
This progress has been brought about chiefly by the application of 
thermal models to the development of low-stress growth configurations. 

A new temperature model was developed which cr.r< analyze the 
thermal effects of much more complex lid and top configurations 

than was possible with the old lumped shield model. 

Growth experiments which supplied input data such as actual 
shield temperature and melt levels were used to verify the modeling 
results. 

Desirable modifications in the melt level-sensing circuitry were 
made in the new experimental web growth furnace, and this furnace has 
been used to carry out growth experiments under steady-state conditions. 

New growth configurations were tested in long growth runs at 
Westinghouse AESD which produced wider, lower stress and higher quality 
web crystals than designs previously used. 
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2. INTRODUCTION 


Silicon dendritic web is a single-crystal silicon ribbon 
material which provides substantial advantages for low-cost manufacture 
of solar cells. A significant feature of the process is the growth from 
a melt of silicon without constraining dies, resulting in an oriented 
single-crystal ribbon having excellent surface features. In common with 
other more classical processes such as Czochralski growth, impurity 
rejection into the melt permits the use cf less pure "solar grade" 
starting material without signif icantly affecting cell performance. A 
unique property of the dendritic web process is the growth of long 
ribbons of controllable width and thickness which not only facilitates 
automation of subsequent processing . nto solar cells, but also results 
in high material utilization since cutting and polishing are not 
required . 

During the previous program (DOE/JPL Contract No. 954654), most 
of the component elements for the reproducible and steady-state growth 
of high-quality web crystals were developed and demonstrated. Area 

n 

throughputs greater than 25 cnr/rain were demonstrated for short periods 
of time. Melt replenishment for periods of up to 17 hours (a one-day 
growth cycle) was demonstrated. Thermal models were developed for 
calculating temperature distributions in the web crystal as a function 
of configuration parameters. 

On the present contract and during the previous reporting 
period, three broad areas of work were emphasized: 

1. The development of thermal stress models in order to 

understand the detailed parameters which generate buckling 
stresses. The model can then be usei to guide the design of 
improved low-stress web growth configurations for 
experimental testing. 
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2. Experiments to increase our understanding of the effects of 
various parameters on the web growth process. 

3. The construction of an experimental web growth machine which 
contains in a single unit all the mechanical and electronic 
features developed previously so that experiments can be 
carried out under tightly controlled conditions. 

Thus, the principal objectives of this work were to expand our 
knowledge and understanding of both the theoretical and experimental 
aspects of the web growth process to provide a solid base for 
substantial improvements in both area throughput and web crystal quality 
and to develop the tools necessary to carry out this objective. 

Although these efforts continued during this reporting period, 
the emphasis shifted to the application of the tools and knowledge which 
had been developed in the previous period to the design and experimental 
verification of low stress web growth configurations. 
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3. PROGRESS IN WEB GROWTH RESEARCH 


3 . 1 Temperature and Stress Modeling 
3*1.1 Introduction 

The overall goal of the present program can be described simply 
as the development of the dendritic web growth process to give greater 
area throughput of silicon rxbwn. One of the most important 
requirements for satisfying this goal is the reduction of thermally 
generated stresses in the web, which translates into creating the proper 
temperature profile along the length of the growing web. Both aspects 
of area throughput, growth velocity and ribbon width, are dependent 
directly or indirectly on the temperature distribution in the ribbon. 

The growth velocity for any specified crystal thickness is 
obviously dependent on the temperature profile, since the, removal of the 
latent heat of fusion is directly related to the temperature gradient at 
the growth front. Although the latent heat can be dissipated through 
both the liquid and solid, our present concern is primarily with the 
heat lost through the solid web; the heat lost to the supercooled melt 
depends on parameters other than those which we will be onsiderlng 
here . 

The width of the web crystal depends in a somewhat more indirect 
fashion on the temperature profile through thermal stresses. These 
stresses are manifested in two effects: residual stress which is 

present in the grown material and buckling of the ribbon which occurs 
curing growth itself. Of these two stress effects, the ’suckling 
phenomenon represents the most severe limitation to the growth of wide 
ribbon since both the magnitude of the stress and the stiffness of the 
ribbon depend strongly on the width of the crystal. Usually a dendritic 
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web crystal is seeded at a narrower width than its desired final 
dimension and then widens as it grows until either steady state is 
achieved (width-limited growth) or until .something occurs which 
terminates growth. Buckling is such an event and is manifested by a 
sudden deformation of the growing ribbon which can, at worst, cause the 
crystal to pull free of the melt and, at best, make the continued growth 
of single crystal difficult. Earlier observations of buckled web have 
shown that material grown immediately prior to the buckling event can be 
essentially free of residual stress, so that residual stress and 
buckling can be ascribed to two somewhat different causes. Once 
buckling has occurred, however, the material usually has a large 
residual stress, although some reoent observations suggest that this may 
not necessarily be the case. It is important to realize that a zero 
stress condition is not required to achieve the goals of the dendritic 
web program. The residual stress observed in ribbon material such as 
dendritic web is the result of thermal stresses exceeding the material’s 
yield point at some time during growth. When the critical yield stress 
is exceeded, plastic deformation occurs through the motion or generation 
of dislocations, resulting in a strained lattice when the material has 
cooled to room temperature. If the critical yield stress is not 
exceeded, then an unstrained crystal can be grown. Although the plastic 
deformation (or, more precisely, the visco-plastic deformation) process 
is extremely complex, nevertheless extremely perfect dendritic web can 
be grown even though the thermal stresses were certainly not zero.^^ 
What is required is to have a thermal stress distribution where the 
critical yield stress (however defined) is not exceeded. 

A similar situation exists with respect to buckling that results 
from thermal stress. In that case, the critical condition depends on 
the ribbon width and thickness (stiffness) as well as on the thermal 
stress itself. As long as the combination of stiffness and stress 
parameters do not exceed some critical condition, flat ribbon can be 
grown. 
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At an earlier point in the development of the dendritic web 
technique, a purely empirical approach was used to design low-stress 
growth configurations. Some progress was made but the rate was very 
slow. Then, compute 1 "’ 3 d e 1 s were developed to calculate the temperature 
profile in a growing web and to calculate the stresses generated by 
these profiles. Most recently, the buckling ^ability of the web was 
calculated, i.e., under what conditions of st:,!ss and stiffness (width, 
^thickness) would the web switch from growth as a flat ribbon to some 


curved, buckled shape 


( 2 ) 


The advantage of this modeling appro^h is that many weeks of 
growth experiments, including hardware fabrication as well as actual 
growth, can be simulated by relatively v Mmplr,, M inexpensive computer 
calculations. The disadvantage, of ^urse, ■lij-.’That the modeling is only 
an approximation of reality and the predictions of a model can be 
trusted only to the extent that the model itself represents the real 
world. In the present case, the modeling has been shown to give a very 
good representation of the web growth results as long as the growth 
geometry can be adequately represented. 


The code which calculates the temperature profile, could be 
verified both directly and indirectly. The direct verification was 
obtained through predictions of the relation between web growth velocity 
and ribbon thickness. The good correlation between the results 
calculated for a real growth system and the measured results for that 
system essentially verified the temperature distribution near the growth 
front. The adequacy of the results for temperatures further from the 
growth front was verified by the prediction of buckling behavior which 
agreed well with the observed web buckling. 

Although the temperature calculation code was successful in 
describing the behavior of relatively simple growth systems — a 
crucible with a lid and one or two radiation shields — it appeared to 
be the weakest element in the analysis for more complex geometries. 

Newer growth geometries had more complicated lid designs and three. 


four, or even more radiation shields. The model was clearly an 
inadequate representation of such a system. The design of a new 
temperature calculation code is described later in this report. 

Before discussing the current year's results in detail, it is 
perhaps well to briefly review the philosophy and approach used for the 
application of the models. A schematic representation of the modeling 
process is shown in Figure 1. In most instances the modeling sequence 
begins with the definition of some growth configuration: the shape and 

temperature of the susceptor lid and the temperature, location, and slot 
size of the radiation shields above the lid. These parameters may 
represent a real growth system being studied in the laboratory or they 
may represent the design of a proposed improved growth system. In the 
former case, real dimensions and measured temperatures are available; 
for proposed designs, obviously arbitrary values can be used. However, 
sint’t toe hardware may eventually be fabricated, physically reasonable 
values are usually chosen. This does not present much of a limitation 
since a wide range of dimensions and element temperatures are attainable 
as determined from past experience. The actual choice of dimensions and 
temperatures for a new proposed design is to some extent empirical in 
that it is based on analytical and experimental results obtained from 
previous designs. 

Once the web temperature profile has been calculated for a 
growth configuration, two possible courses of action are possible. 

First, an iteration on the design could be made on the basis of the 
temperature characteristics (the (aT)" loop shown in Figure 1). 
Alternatively, the temperature data can be used as input to the finite 
element code that calculates thermal stresses in the ribbon, and 
iteration on the proposed geometry can be done on the basis of the 
calculated stress distribution. Alternatively, the stress distribution 
can itself be used as an input to the modeling code that calculates the 
buckling eigenvalues for the web. 
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Not all the temperature profiles used as input to the stress 
model resulted from analysis of web growth systems. Some temperature 
distributions were defined from purely mathematical considerations to 
illustrate isolated features of the "real" profiles. For example, a 
temperature profile which should produce a constant stress in a long 
ribbon was used to illustrate the effect of the free boundary at the 
growth front. The results of several of these "synthetic" profile, runs 
will be discussed later. 

In the following sections, we shall discuss the development of 
an improved model for calculating web temperature profiles, the results 
of stress calculations for several "synthetic" temperature profiles, and 
finally the application of the models to the development of improved 
growth configurations. 

3.1.2 Model Development 
3. 1.2.1 Introduction 

The modeling of the thermal stress and buckling of dendritic web 
crystals in fact involves three distinct computer codes: 1) a model 

which calculates the temperature distribution along the growing ribbon; 
2) a model which uses the temperature distribution as input and then 
calculates the thermal stress; and 3) a model which uses the thermal 
stress distribution as input and then calculates the stability of the 
ribbon with respect to buckling. With such a sequential type of 
analysis, it is of obvious importance to have as good as possible a 
model in the first stage. In fact, the second and third stages of 
modeling utilize a general, proprietary finite element code, WECAN, 
which has been shown to be more than sufficiently accurate for our 
present purposes. The present work concerns improvement of the 
calculation of the temperature distribution in the web which forms the 
basis for the subsequent calculations. 

The model which calculates the temperature distribution along 
the growing silicon ribbon has two distinct routines: a routine which 
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calculates the radiative heat exchange between a point on the ribbon and 
its environment, and a second routine that integrates the heat conduc- 
tion equation. The problems which occasioned the present work lay not 
in the integration routine, which had been upgraded earlier, but in the 
geometrical factors. 

For several years, the geometry shown in Figure 2 has been used 
to evaluate the radiation interchange with the ribbon. This model quite 
adequately represented the simpler growth configurations used until 
recently, which had only one or two shields above the lid. Currently, 
growth configurations with multiple, and occasionally widely spaced, 
shields are being used, and representation of such a configuration by an 
isothermal rectangle is suspect. During the reporting period, we have 
developed a new computer code which allows calculation of the radiative 
interchange between the web and a number of shields, spacers, etc. so 
that we have greatly improved the geometric resolution of the model. 

As in the previous model, we have neglected other heat transfer 
mechanisms such as gaseous conduction and convection. Furthermore, we 
still consider that the absorbance of the silicon ribbon is equal to the 
emittance, and that the shields, lids, etc. are black bodies so that 
multiple reflective transfer is neglected. While such refinements might 
be desired from a theoretical standpoint, the success of the previous 
modeling suggests that they are second order effects. 

3. 1,2. 2 Review of Radiative Heat Transfer 

Consider the radiative transfer between a single shield at 
Temperature T^ and an element dx of the web at the point x at 
temperature T(x) (which is to be determined) . Figure 3 illustrates the 
shield subtending an angle from 0^ to 0 i2 f rom the normal to the 
element dx. The radiative flux density from the shield to the point x 
on the web is 

’ ^ , 12 

2 ° T i V 036 * (1) 

il 
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where o is the Stefan-Boltzraan constant. The radiative flux density 
from the point x to the shield is 


1 A r 

T eo T(x) / 


12 


cosO d0 


il 


( 2 ) 


where e „.s the emissivity of silicon web. Since we are assuming the 
absorbance of the web equals its emittance, the new flux density 
transferring from x is 


q = ea [T(x) 4 - T. 4 ](sin0 i2 - sin0 u )/2 


(3) 


The net flux density transfer from any number of elements is the 
sum of terms like che one above. The essential fact is that once the 
angular limits and temperature of each element are determined, q may be 
simply evaluated for substitution in the heat transfer Equation A: 


P Cu 


dT_ 

dx 


dx'T dx"' 


_ ii 

t 


(A) 


where p = density; C = specific heat; u = pull velocity; a = 318 W/cm, 
the coefficient of the temperature dependence of thermal conductivity; 
and t = web thickness. As before, this equation may be integrated to 
determine the temperature distribution T(x) along the web. We need only 
evaluate the geometric factors of Equation 3. 


3. 1.2. 3 New Geometry 

We examine now the original geometry of Figure 1 to see how it 
can be systematized for extending its number of elements. A possible 
systematic arrangement is illustrated in Figure A; by numbering the 
ambient regions in sequence with the lid and shields, they may be 
treated in the same fashion as the other elements. The melt region and 
the crucible wall are assumed to be at the same temperature and are 
assigned the number "0”. In this case, the ambient region above the top 
shield would be numbered element 5. The lid is divided into two 
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elements because of the bevel. Since the heat transfer depends only on 
angular subtense (reflections are neglected), the beveled lid is 
equivalent to a stepped lid as indicated by the dotted lines. In the 
previous geometry, the whole lid was isothermal, but now, if desired, 
different temperatures can be assigned to the top and bottom elements of 
the lid. In any case, replacing the bevel with a step makes all the 
elements rectangular in shape. Since the geometries of interest have 
the back edges of the elements lining up, the ambient spaces can be 
considered to be comprised of rectangles such as element 3. Thus, the 
geometry is completely determined by the input of the upper corner of 
each element indicated by the large, black points. 

The four numbered points on the web are the boundaries of 
regions with different views; for example, for points between 3 and 4, 
all elements are in view, while above 4, the ambient region 3 i~ 
invisible. These points along the x-axis can be systematically 
determined in the following manner. First, the lower "shade" of each 
element is found by drawing lines through the inner upper corner of each 
element below it. The highest intersection of these lines on the x-axis 
is the lower shade. For example, the point 1 is the lower shade of the 
4th element; element 4 is invisible for points on the web below point 1. 

After all the lower shades are found (some elements may not have 
any), then the upper shades are found. The lower inner corner of an 
element is connected to the lower inner corner of each element above 
it. The lowest intersection point of each of these lines with the 
x-axis is the upper shade of the element in question. For example, 
element 2 does not have an upper shade, while element 3's upper shade is 
point 4 on the web; element 3 is invisible above point 4 on the web. 

After all the upper and lower shades are found, these points are 
ordered on the x-axis. The region between two cuch points is assigned a 
vector, the ith component of which is either one or sero depending on 
whether or not the ith element is visible to points la this region. 



Figure 5 illustrates a possible geometric configuration for the 
new model. Numerical integration of the heat conduction equation 
proceeds up the silicon web from the melt surface. By knowing in which 
region between upper and lower shades the point of integration lies, we 
know which elements are visible to it. For the point illustrated, 
elements 0 (the melt), 4, 6, 7, 8, and 15 (ambient) are visible. There 
is one additional piece of information needed before the integration can 
proceed: the element opposite the point of integration, element 7 in 

this case. For the visible elements numbered 7 or less, the lower inner 
corners are used to limit the viewing angles, while for elements 7 or 
greater the upper inner corners are used as illustrated. The viewing 
regions are labeled 1 through 6; these determine the proper coefficient 
q of the heat conduction equation during integration. 

Once the geometrical heat transfer coefficients have been calcu- 
lated, Equation 4 can be numerically integrated to yield the temperature 
distribution along the length of the web crystal. A more detailed 
description of the computer code is presented as an appendix to this 
report. 

3. 1.2. 4 Verification of the New Model 

One of the modeling runs of the J419 configuration, Case 9-11C, 
was used to verify the new radiation transfer model. The first test was 
to run the identical configuration, i.e., a lumped shield model as in 
Figure 2, to verify the operation of the routine which evaluates the 
shade points, etc. The results were excellent; the new program gave the 
identical temperature profile as the older program. A second geometry 
was then evaluated in which both the lid and the shield "block" were 
beveled. This is a relatively small change in the geometry, but one 
that could not be accommodated by the older model. It would be expected 
that the results be similar, but not identical, and that was indeed the 
case. 
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Figure 5. Viewing regions from a point on the web. 
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Figures 6 and 7 illustrate the two cases.* However, instead of 
presenting the temperature profiles themselves, the figures show the 
second derivative of temperature times thermal expansion coefficient, 
i.e., (aT)". This parameter was chosen because it accentuates the 
curvature of the temperature which can basically be considered as the 
cause of thermal stress. Examination of the two figures shows that 
indeed the two profiles are almost the same, although there are some 
slight differences due to the more complicated profile. 

3.1.3 Operation cf the Model 

Input options that are available for the model fall into two 
categories: selection of the operating options of the model such as 

type of output, number of cases, mode of integration, etc. and, second, 
the parameters of the configuration being analyzed. A description of 
the individual data cards as well as a complete listing of the program 
are included as an appendix to this report. 

Output options available for the program fall into three 
categories: 1) numerical output, 2) graphics output, and 3) punched 

card output. The first and third output types were previously 
available, the second output was added when it became apparent that the 
second derivative of the thermal expansion coefficient times the 
temperature was a useful parameter for assessing the effect of changes 
in the shield elements on the thermal stress. 

An example of the numerical output is stiown in Figure 8. The 
first column lists the position on the web (the growth front is at 0 
cm). The second column gives the calculated temperature. The third and 
fourth columns are the first and second derivatives of the temperature, 
and the fifth column is (aT)". All distances are in centimeters and 


*These figures were generated by a computer plotting routine which has 
been added to the temperature distribution model. Not only is (aT)" 
plotted, but also the lid and shield geometry is represented (heavy 
lines in the figures). 
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Figure 7. Model geometry and (oc T ) ” for run 9-11CA (beveled elements). 
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temperatures are in °K. In addition to this data, two other sets of 
data are listed at the top of the printout. The first, labeled VV, is 
the partial growth velocity of the web resulting from the heat lost from 
the web crystal itself. This differs from the total growth velocity 
which is observed in experiments hv a contribution which results from 
some of the latent heat being dissipated to the supercooled melt.(^) An 
estimate of this total growth velocity is included as one of the input 
parameters to the model. 

Tne second set of data included is an estimate of the critical 
yield stress for the web corresponding to the first five temperature 
points. The yield stress is simply an estimate based on an empirical 
equation from the work of Graham et 

u Y p = 2.57 x 10" u exp(49459/T) Mdyn/cm 2 (5) 

These critical-yield stress values can be used later when actual stress 
distributions are generated by the finite element WECAN calculations. 
Although the viscoelastic phenomena responsible for the observed 
residual stress in the crystal*', is extremely complex, these data provide 
a first estimate of possible effects. The first four columns in the 
balance of the numerical data are more or less self-explanatory in being 
the position, the temperature, and its first two derivatives. The last 
column, (aT)", is of importance since it is basically the generating 
function for the thermal stress. Although the relationship between 
(aT)" and stress is not straightforward near the interface where the 
free boundary exerts a strong influence, it is relatively simple 
elsewhere where the Boley and Weiner approximation^^ applies: 

°x = 6* < ‘ 3y2 “ ^ aT ^" 

where E is Young’s modulus and w is the half width of the ribbon. 

Through the relationship of T" with the heat loss from the ribbon, it is 
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possible to correlate the behavior of this parameter with changes in the 
geometry of the lids and shields. 

Realization of the usefulness of the (aT) M parameter led to the 
addition of a graphics output capability to the code. Since the end 
application ” oc ’ to compare the (aT)" function with the geometry being 
analyzed, the graphics output presents both a representation of the lid 
and shields and the concomitant (aT)” curve. If the geometry does not 
change (lid, shield, and interface position), then up to three different 
(aT)" curves can be represented on the same plot. A representative 
output is shown in Figures 6 and 7; the plot of the lid and shields has 
been accented for the sake of emphasis. 

The final output option from the program is a set of punched 
cards giving the nodal temperatures for use with the WECAN stress 
calculations. In this option, there are several sub-options depending 
on whether a two-dimensional stress calculation or a three-dimensional 
buckling calculation is to be done. Further, in the two-dimensional 
calculation, there is a choice as to whether quadratic or cubic elements 
are employed. Also, some choice is available as to some of the 
geometric features of the finite element grid. 

Thus, the program to calculate the web temperature distribution 
has a great deal of flexibi Lity in both the lid and shield configuration 
which can be analyzed, and in the options for presenting the results of 
its calculations. The flexibility of both the input and output greatly 
enhance the ability of the program to assist in the analysis and design 
of dendritic web growth configurations. 

3.1. A Stress Model Development 

Although most of the development of the models themselves was 
concerned with the temperature code, some investigations of the stress 
and buckling models were also made. The concern was not so much with 
the computer code itself, which is well validated, as with the 
application of the proper boundary conditions, especially with the 
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buckling calculations. The question gained importance as wider and 
wider crystals were being grown and the growth systems analyzed. 
Calculations with "synthetic" temperature profiles suggest that the "end 
effect" due to the free boundary at the end of the web extends longitu- 
dinally for a distance about equal to the full width of the ribbon. 

Although this should be a good approximation to reality for the melt end 
of the crystal, it could lead to an underestimate of buckling when also 
used at the "ribbon end" of the finite element mesh, especially when the 
crystal width approached 5 cm and the mesh length was only 10 cm. A 
longer mesh could of course be used, but to make much difference the 
increased length would be significant and would lead to increased 
calculation time. A simple change in the boundary conditions could 
possibly accomplish the same end and would not affect computer usage. 

In calculating the stresses resulting from a. given temperature 
profile, appropriate temperatures are assigned to the nodes of a mesh of 
elements. In addition, boundary conditions must be specified for the 
edges of the mesh. At the boundary representing the growth front, the 
appropriate condition is the free boundary a x = 0, i.e., no longitudinal 
stress. On the free edge of the strip, the condition is o y = 0, i.e., 
no lateral stress. The center of the strip is a line of symmetry so 
that the condition there is zero y displacement u^ = 0. The uncertainty 
arises at the end of the mesh which "connects" to the rest of the 
ribbon; previous practice has been to consider this to be another 
traction-free boundary with c = 0 as at the melt end of the mesh. A 
possible alternative condition would be to require that all x-displaceraents 
be the same at the boundary, i.e., u x = const. 

The two possibilities were tested by calculating the buckling 
eigenvalues for a model of the J460 configuration, which was known to be 
capable of growing undeformed crystals over 40 mm wide at thicknesses of 
150 pm or less. A 39.5 mm width was assumed for the model (37 mm web 
with 1.25 mm dendrites on each edge). The free-boundary condition 
predicted a cl Ideal width of 39.9 mm for buckling, which is slightly 
narrower than actually observed, while the constant displacement 
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condition predicted a critical width of 35.1 mm. Obviously, the free 
boundary condition is much closer to observation and in fact is slightly 
conservative. Examination of the stress distributions showed that the 
"free-boundary" condition was probably much closer to the probable 
actual stress than the "constant-displacement" condition, in agreement 
with the buckling predictions. It would appear therefore, that our 
present application of the stress and buckling codes is adequate for 
evaluating material up to about 45 or 50 mm; but when wider growth is 
predicted, then the finite element mesh must be extended to preserve 
good agreement between the model and reality. 

Another topic which was investigated was the relation between 
ribbon width and stress magnitude. Usually a 27 mm crystal is assumed 
for calculation of stress when evaluating a temperature profile with a 
two-dimensional stress model. This width has been appropriate to real 
crystal dimensions and by holding the width constant, comparison of 
various temperature profiles was facilitated. With the growth of wider 
and wider crystals, however, we wanted to know how stress might vary 
with width. 

A temperature profile was calculated for one of the J460 type 

configurations and two-dimensional stress calculations were made for 

17 mm, 27 mm, and 39.5 mm models. Surprisingly, the y-stress 

distribution near the growth front varied only slightly with width, 

especially for the two wider models. The x-stresses further from the 

growth front varied more strongly with width, but the relationship was 
13 2 

more nearly W ' than W as expected from simple theory. 

These results are consistent with some observations of web 
growth in the laboratory. When a poor growth configuration is used, 
even narrow crystals will have relatively large residual stress, which 
would evidently be expected if the y-stress near the growth front is the 
causitive factor. 
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3.1.5 Synthetic Temperature Profiles 


Previous studies of thermal stress and buckling have used 
mathematical models to evaJuate the performance of real, or at least 
proposed, lid and shield configurations. This type of analysis has been 
continued during the present period but, in addition, stresses were 
computed for several temperature distributions which were synthesized 
purely on the basis of their mathematical properties. These "synthetic" 
temperature profiles were evaluated to illustrate the requirements of a 
low-stress profile. 

Boley and Wiener^) have published an often-cited series 
solution for the stress components in a uniform ribbon having a 
temperature distribution varying only with x (the coordinate along the 
length of the ribbon). The first term of their result for o xx is 


a 

xx 



2 

r 


[6a] 


where a is the thermal expansivity (assumed to be constant), E is 
Young's modulus (also constant), w is the half width of the ribbon, and 
y is the width coordinate. The prediction of equation 6a, that a xx = 0 
if T" = 0, does not hold if the thermoelastic properties (a or E) vary 
with temperature; for our calculations we assume that E is constant but 
allow a to depend on temperature. The condition, for zero stress is then 
(aT)" = 0, and for constant stress, Equation 6a becomes 

(aT)” = 6o xx /E(3w 2 - y 2 ) [6b ] 


By assuming an analytic form for a(T), e.g., a = a Q + a^T, Equation 6b 
can be solved to identify temperature distributions giving constant 
stress in the ribbon. A series of such temperature distributions is 
shown in Figure 9, where the parameter labeling the curves is the 
difference in stress between the ribbon center and ribbon edge. 
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| Figure 9. Constant stress temperature profiles. 
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The series solution from which Equation 6a was derived only 
applies when end effects are negligible and when T" (or(aT)") is not 
changing rapidly. Foth of these conditions are violated near the growth 
front of the dendri';ic web. Since essentially zero tension is applied 
to the growing interface, the growth front is a traction-free 
boundary. Further, because of the rapidly changing high temperature of 
the ribbon in this region, T" also changes rapidly. In order to develop 
at least a sami-quantitative understanding of the effect of these two 
conditions on the thermal stresses, two "synthetic” temperature profiles 
were used as input data for the WECAN stress calculations. 

To illustrate the effect of the traction-free boundary, a 
temperature profile was used that would generate a constant stress 
according to Equation 2. The traction-free boundary condition was 
intrinsic to the WECAN analysis used for the calculation. 

The modeling of this profile was actually run twice. The first 
run gave results which were very close to the predictions of Equation 6b 
but differed in some small details. It was then realized that the mesh 
used in the stress computation included elements which represented the 
additional stiffening due to the bounding dendrites. The same 
temperature profile was then rerun using a mesh representing a uniform 
ribbon without dendrites, and the results were in exact agreement with 
the theory. 

The result is shown in Figure 10, which depicts A a along the 
ribbon. This stress parameter was chosen since it reduces the "grain" 
inherent in the finite element calculations. From simple inspection of 
the curve, it can be seen that the end effect of the free boundary 
extends for a distance approximately equal to the ribbon width (2.7 cm 
in this calculation). Although the fit is not exact, the effect of the 
boundary seems to be approximately a complementary Gaussian with a 
characteristic length equal to the ribbon width. This result shows that 
the details of the (otT)" curve near the growth front are mitigated to a 
large extent by the end effect. This may be of importance in future 
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lid/shield design efforts in that the design itself may be separable 
into two regions: one near the growth front and one further up the 

crystal . 

The second synthetic temperature profile was intended to show 
how stresses caused by changes in (aT)” affect adjacent regions of the 
crystal. The temperature profile used for the analysis was composed of 
two zero stress segments with the transition at the midpoint of the 
finite element mesh. Thus, (otT)’’ = 0 everywhere, but all the higher 
derivatives were essentially infinite at the midpoint to give an 
effective delta function for thermal stress generation. 

The results of the WECAN calculations are shown in Figure 11; 

11a is the temperature distribution used as input for the WECAN 
analysis, and lib is the resulting o x and 0 ^. The effect of 
discontinuity in the curvature of the profile obviously spreads over a 
significant region of the crystal. Again, the o x distribution is 
approximately a Gaussian with a characteristic length of the ribbon 
half-width. 

The importance of the results in this run lie in understanding 
the interaction of different regions of the (aT)" curve. In the 
previous analysis of web growth configurations, it was obvious that 
there was a connection between the stress components and the (aT)" 
values, but it was certainly not as simple as implied by Equation 6a. 

The present result indicates that an averaging (or perhaps a 
convolution) analysis is more appropriate but further gives some 
indication as to the magnitude of the characteristic distance. 

The third synthetic temperature profile assumed an exponential 
decay of the (aT)" profile from its value near the interface. This 
particular function was chosen since it is a reasonable approximation to 
the behavior observed in models of any realistic lid/shield 
configurations, although in those cases the exponential behavior is 
found only in the first 5 mm or so. The stress distribution from the 
synthetic temperature case was in reasonable quantitative agreement with 
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temperature profiles. 
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the stress distributions in the interface region of the more realistic 
cases, although of course the more distant features were quite dif- 
ferent. In both the synthetic and realistic cases, the actual x-stress 
magnitude near the boundary was far smaller than might be anticipated, 
in agreement with the boundary effect of the constant stress case. 

These results indicate that the initial (aT)" peak is responsible for 
much of the thermal stress behavior at the interface, but also that 
there are other factors such as the free-boundary effect which must be 
considered. It is also not yet clear to what effect the peak height and 
the characteristic decay length of the exponential are involved in the 
stress generation. These factors need some additional runs to clarify 
the behavior. The final results should be applicable to the design of 
lid configurations for faster growth with lower residual stress. 

3.1.6 Development of New Configurations 

Although the model development and the studies of synthetic 
temperature profiles described in the preceding sections have some 
intrinsic importance, their primary purpose is to provide the tools and 
guidance for the development of improved dendritic web growth systems. 
During the reporting period, two principal designs were developed: the 

J419 and J460 series of growth hardware. Each proved to be a distinct 
improvement over preceding designs and each has gone through a series of 
design modifications for use in melt replenished, steady-state growth. 

J419. As indicated by the design designation (run J419 was the 
first experimental evaluation of the configuration), J419 was the first 
of the two designs. Development of the design followed from the 
hypothesis that buckling in dendritic web growth was related to a large 
peak in the x-stress distribution which occurred several centimeters 
above the growth front. This peak in the stress distribution in turn 
appeared to be related to a small maximum in the (aT)" distribution 
calculated by the temperature model. These features are illustrated in 
Figures 12 and 13 for a temperature profile representative of the J98M3 
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configuration, one of the best prior configurations.^) The maximum Aa x 

l» o x (center - o x (edge)] is 615 Mdyn/cm^ and the associated (aT)" is 

—3 —2 

1.2 x 10 J cm . The goal of the development work was to reduce these 
parameters . 

The design of a new configuration is based on prior experience, 
both analytical and experimental, by observing the relationship between 
growth performance and the characteristics of the configurations as 
calculated from the models. Experience with various modeling runs, both 
for real growth configurations and with the "synthetic" temperature 
profiles, suggested that some additional localized heating was needed 
near the growth front while maintaining the loss to furnace ambient. In 
addition, losses further along the crystal should be reduced to lower 
the large x-stress peak. These design requirements translated into a 
relatively thick, hot lid with a relatively tall shield stack having 
shield slots of increasing width. 

The resulting configuration is shown in Figure 14. The lid 

itself is thicker — 0.5 inch instead of the .38 inch or .25 inch lid 

used in previous configurations. Surmounting the lid are four shields 

with slot openings arranged so that the web at the growth front has a 

25.6° clear view of the cold furnace interior. This configuration was 

analyzed for thermal stresses using lid and shield temperatures deduced 

from measurements on other configurations. The temperature calculations 

indicated a "far peak" value for (aT)" of about 8.1 x 10"^ cm"^ and the 

Aa x profile, shown in Figure 15, shows that the accompanying stress has 

o 

been reduced to about 360 Mdyn/cm . 

On the basis of the calculated thermal stresses, a full buckling 
analysis was performed on the model. The results Indicated that a 150 
ym thick web crystal would not buckle at widths narrower than 38 mm. 

This was a substantial improvement over previous designs and the 
decision was made to fabricate a set of lids and shields that would give 
a real growth system representing the design analyzed with the models. 
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Figure 12. Temperature model of J98M3 configuration. 
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' Figure 13. Delta x-stress for J98M3 configuration. 
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Figure 15. Delta x-stress for J419 configuration. 
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The initial run of the configuration included thermocouples in 
the lid, bottom shield, and top shield. The measured temperauui.es were 
consistent with the temperatures assumed for the stress and buckling. 
More important, however, was the result that in the initial run an 
unbuckled crystal was grown to a width of 35 mm, although the thickness 
was only 125 yra. In subsequent runs, thicker web (190 ym) was grown to 
a width of 44 mm without buckling and to 49.9 ram with only slight 
deformation. These width and thickness values are in almost perfect 
agreement with the buckling predictions. 

J460. Although the J419 was a marked improvement over previous 
configurations, it was felt that further improvement could be made 
following the same design philosophy. A thicker lid would hopefully 
improve the temperature distribution near the growth front to give low 
residual stress, while an even higher shield stack would reduce the 
stress peak away from the interface to reduce the "buckling" stress. In 
choosing a geometry to model, a recessed lid of the J352 type^^ was 
used for LI with a solid lid for L2; the recessed lid would minimize gas 
conductivity effects while not changing the radiative heat transfer. 
Thus, a real system would more closely approach the assumptions in the 
modeling. 

The results of the temperature profile calculations are shown as 
a (aT)" plot in Figure 16. The recessed lid geometry is apparent from 
the geometry plot (accented) which is part of the figure. In this 

/ O 

calculation, the (aT)" maximum has been reduced to 2.3 x 10” cm and 
the resulting Ao x is only 195 Mdyn/cm / as seen in Figure 17. The growth 
experience with this configuration confirms the low-stress aspects 
predicted by the models. 

In summary, the use of temperature and stress models to identify 
the sources of thermal stress in dendritic web growth systems has led to 

marked improvement in the growth hardware. Very early configurations 

2 

had stress peaks as large as 1200 to 1500 Mdyn/cm , the J98M3 con- 

O 

figuration had a peak of about 615 Mdyn/cm ; the J419 peak was about 
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Figure 16. Temperature profile results for J460 configuration. 
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350 Mdyn/cm , and in the J460 configuration the peak, was only 195 

O 

Mdyn/cm . With each reduction in peak height, wider web was grown. 

Even further improvement will be possible. 

3.2 Operation of the New Experimental Web Growth Facility 

During the previous reporting period a new web growth facility, 
the N-furnace (shown in Figure 18), was constructed . This facility 
includes all the features required for sustained steady-state web 
growth. It permits operation under conditions of constant melt level, 
constant temperature, constant crystal width, thickness, and speed of 
growth and has provisions for programmed start of web growth. 

During the early part of this reporting period, growth runs in 
the N-furnace were dedicated to achieving the proper furnace adjustments 
for continuous feeding of a standard width-limiting J98M3A lid and 
shield assembly combined with the elongated crucible. This configur- 
ation, which had been well characterized in the WA and WB furnace 
facilities, although not with controlled melt level, was adapted for 
melt replenishment by the addition of feed and laser holes. 

When initial testing began, installation of the melt level- 
sensing circuitry had not yet been completed, so that replenishment was 
carried out with manually set feed rates. The first step was to 
establish the end shield settings for various replenishment rates so 
that sufficiently high temperatures were maintained in the feed 
compartment of the crucible to melt the pellets at the required rate. 

At the same time, growth behavior was tested at these shield settings in 
order to establish what feed rates could be achieved without 
compromising web growth, i.e., the amount of shielding that could be 
used without seriously perturbing the temperature distribution in the 
growth region of the melt. These experiments indicated that for web 
widths of about 2.5 cm, full replenishment rates could be accomplished 
without difficulty, but that the shielding required for feed rates 
needed for 3 cm or greater width could begin to affect the temperature 
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Figure 18. New experimental web growth furnace. 
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distribution in the growth region adversely. This problem was corrected 
by a minor modification in the lid geometry. 

During this period, installation of the melt level-sensing 
circuitry was completed and tested during growth experiments. Although 
the melt level control system components all functioned, it became 
apparent that the system as a whole was not properly tuned to the growth 
requirements. Several changes were made in the control loop which 
brought the response time and loop gain into better agreement with the 
control requirements, and excellent results were then obtained with 
actual web growth. In one run, several crystals were grown under 
automatic level control. During the whole course of the run, all the 
crystals maintained a constant thickness within about 4 pm (at one 
constant growth speed) and even the hold temperature did not change by 
more than a few microvolts. The extreme stability of the growth 
conditions would indicate that the level control system was operating 
even better than required. 

The modifications made to the N-furnace melt level control 
system to better meet growth requirements were three in number. First, 
the time constant of the control loop was made adjustable so that 
integration periods of 1, 20, 100, and 400 seconds could be used. The 
faster response times are used for the initial set up of the equipment, 
buL the longer times are needed to avoid unnecessary fluctuations in the 
final feed rates. With the longer integration times, a higher loop gain 
can also be used and the response gain of the sensing circuit was 
increased by a factor of about ten. As a result, the control 
capabilities of the system are easily of the order of 100 pm or better 
in melt height sensing. Finally, a selectable limiting speed was placed 
on the pellet feed motor in order to prevent overfeeding during 
transient conditions. 

At this point, the N-furnace became a web growth research tool 
capable of generating growth conditions controlled to a degree not 
previously available. 


43 


3.3 Experimental Web Growth 

3.3.1 Introduction 

The process of developing a functional new growth configuration 
follows a three-stage progression. The first stage is the thermal 
modeling, the results of which generate a design for a growth 
configuration. Stress and buckling models are applied as appropriate 
when the (aT)" results warrant the effort. The design is then 
fabricated into hardware and tested experimentally. The objective of 
this second phase of the progression is to experimentally verify the 
stress behavior predicted by the model, i.e., to find how wide the 
crystal can grow before deformation occurs. At this stage, long slots 
which provide melt temperature profiles compatible with crystal widths 
of about 6.5 cm are used. Lid and shield temperature measurements are 
made to verify consistency with the model, and modifications are made as 
required to product the desired temperature distribution in the vertical 
direction. This stage involves a cross interaction between modeling and 
experiment. The third stage of development is the adaptation of the 
low-stress configuration for semi-automated growth. This phase is 
largely empirical, guided by experience. It involves the incorporation 
of melt replenishment and width control provisions into the design, and 
the determination of optimum shielding and coil position for good growth 
at constant width and melt level, i.e., steady-state growth. Only the 
designs which continue to show promise through the first and second 
phases reach the third phase. 

3.3.2 The J419 and J435 Configurations 

The first lid and top shield configuration in which hardware 
design was generated by the thermal stress modeling results was the 
J419, discussed in Section 3.1 of this report. The stress 
characteristics of crystals grown from this configuration agreed 
extremely well with the predictions of the models, with buckling 
stresses much lower than those obtained with previous configurations. 
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Encouraged by the success of the baseline J419 design, a series 
of variations of this configuration were tested in order to test the 
sensitivity of growth behavior to variations such as top shield spacing 
and minor variations in slot geometry and coil position, etc. Lid and 
shield temperature measurements were made to complement the modeling 
work. Growth parameters evaluated were growth velocity and residual 
stress, in addition to the width at which buckling was initiated. 

In one series of experiments, the J419 was modified to a more 
open geometry by beveling the top lid and using a different shield stack 
order. Although the growth speed was enhanced, there was also a change 
in the character of the residual stress. Whereas the standard J419 
configuration produced crystals with negative residual stress, the 
modified version gave material with positive residual stress. In both 
cases, the magnitude could be either very small or moderate, depending 
on other growth parameters. The higher speed modifications also would 
grow moderately wide material, up to 40 mm, but tended to deform at 
slightly narrower widths than the unmodified J419. 

The earlier modeling results suggested that the residual stress 
in some of the higher speed versions of the configuration could be 
reduced by reducing the heat loss from the web in the interface 
region. One feature of a lid design which would accomplish this would 
be a bevel on the bottom edge of the slot in the lid. Thermally, this 
would allow the web to "see" the hot crucible cavity for a greater 
portion of its length so that the heat loss near the growth front would 
be decreased. 

Four runs (J439-J442) were made in which the top lid was either 
unbeveled or had a 3 mm bevel; the bottom lid had either a 1.5 mm or a 3 
mm bevel. It was found that the bevel on the bottom of the lid had more 
effect on speed than on residual stress, while the bevel on the top lid 
affected the stress more than the speed. 

In general, the only modification of the basic design which had 
a positive effect on growth was a small bevel on the lid slot. Other 
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variations either had little observable effect on growth behavior or 
generated a negative result. For example, opening the slots in the top 
shields increased the growth velocity but caused the crystal to 
degenerate at narrower widths. The model indicated that increasing the 
height of the top shield stack should further reduce buckling stresses, 
but this did not seem to have a significant effect with the J419 lid 
configuration. 

While modifications to the J419 configuration were being tested 
in the J-furnace, experiments in the WA-furnace were directed to testing 
the compatibility of the J419 configuration with the elongated crucible 
which would be required for growth at constant melt level. Adjustable 
end shields were installed. Melt probe data showed that the lateral 
coil position and the end shield height affect the melt profile 
differently and that adjustment of both parameters is necessary to 
obtain a flat, symmetric temperature distribution in the melt. However, 
when the proper melt temperature profile was achieved, it was 
established that the low-stress characteristics of the J419 configura- 
tion were maintained with the elongated crucible configuration. 

At this point it seemed appropriate to combine the wide-growth, 
low-stress capabilities of the J419 configuration with the width- 
limiting behavior of the J98M3A configuration for demonstrating extended 
growth runs. Lids and shields were fabricated which combined these 
features in a hybrid configuration now known as the J435. Four runs 
were made with this configuration in the J-Furnace to evaluate the 
width-limiting capabilities of the design, and the results confirmed a 
steady growth width between 28 mm and 33 mm, depending on the furnace 
parameters . 

Feed holes and laser holes were then added to the J435 lids and 
shields and the hardware installed in the N-Furnace with an elongated 
crucible for melt replenishment. Preliminary experiments focused on 
optimizing the end shield adjustment and the work coll position so that 
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the melt profile in compatible with both steady-state width and melt 
replenishment. 

Several growth runs were made both with and without melt 
replenishment. In one run, a 4.8 meter crystal was grown with the width 
held at 3.1 to 3.3 cm for 3.4 meters. This crystal was grown without 
replenishment due to clogging of the feed hole with oxide, a problem 
later resolved. Several 3 to 4 meter long crystals were grown with 
replenishment, with thickness constant to within a few microns and width 
limited in the 3.1 to 3.3 cm range. 

The J435 configuration was installed in furnaces at the 
Westinghouse Advanced Energy Systems Division (AESD), replacing the 
J98M3A configuration, where it could be tested in long (100 hr) growth 
runs. Long undeformed crystals were produced regularly at widths of 

3.3 cm. Additionaly, these crystals had much lower residual stress than 
crystals grown from the J98M3A configuration, resulting in a significant 
reduction in breakage during cell processing. While the experiments 
».fith the J419 and variations were being carried out, the new temperature 
model was developed. This model can accommodate complex lid slot 
designs and any discrete number or spacing of top shields, enormously 
expanding the range of designs which can be accurately evaluated. 

On the basis of results using the new model, combined with 
experimental correlations and experience with the J419, it was concluded 
that a new lid slot geometry combined with an extended shield stack 
should produce significantly lower buckling stresses than the J419 
configuration. This new design has been designated the J460 
configuration. 

3.3.3 The J460 Configuration 

The J460 configuration has a combination lid composed of a 9.5 mm 
thick J352 lid with a recessed slot for the bottom element (LI) and the 

6.3 mm thick top section (L2) from the J419 configuration. The intent 
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of the recessed section of the lid was to reduce gas conduction effects 
and so to emphasize the radiative transfer assumed in the modeling. 

During the initial experimental evaluation of this design, web 
crystals up to 48 mm in width were grown without obvious evidence of 
buckling. Furthermore, ribbon split measurements indicated very low 
residual stress levels, even on 41 mm wide samples. The latter data 
were very encouraging since the residual stress might be expected to 
increase as the third or fourth power of the ribbon width and thus be 
fairly high in so wide a crystal. 

On the basis of this initial success, the J460 configuration was 
installed in a furnace at the Westinghouse Advanced Energy Systems 
Division (AESD) in order to accelerate the accumulation of growth 
experience. A number of undeformed crystals 5.0 to 5.4 cm in width were 
grown. These crystals were generally terminated not because of 
degradation in quality, but because of melt depletion (continuous 
replenishment was not in use at AESD at that time). Furthermore, the 
high quality of material grown from the J460 is evidenced by the fact 
that solar ceils fabricated from it have consistently shown efficiencies 
in the 15 to 16% range. 

Having evaluated the basic growth behavior of the J460 
configuration, it was then adapted for melt replenishment in order to 
test growth at known melt levels. The top shields were also 
instrumented to measure actual shield temperatures at different coil 
heights for incorporation into the models. A device was constructed to 
accurately measure melt position so that growth parameters such as the 
velocity-thickness relationship could be determined at known growth 
conditions. This work is still continuing. 

The next step was to design a width-limiting version of the 
J460. This design has been completed and will be tested experimentally 
during the next reporting period. 
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3.3.4 Long Growth Runs and Oxide Control 

As part of a parallel web growth development program funded by 
Westinghouse AESD, four runs were made in which two furnaces, one of 
them the N-furnace, were operated 24 hours per day for up to 103 hour 
long periods to evaluate possible practical problems with continuous 
melt replenishment for extended times. No fundamental difficulties were 
encountered with the operation of the equipment. However, some problems 
with oxide accumulation in the melt sensor laser holes and the feed 
holes in the top shields were identified. In the case of the feed 
holes, the dropping pellets served to keep the feed hole clear, but if 
replenishment was interrupted, oxide tended to accumulate until the hole 
was blocked. 

The oxide accumulation in the laser holes was eliminated by 
fairly straightforward changes in the geometry of the holes in the 
shield stack. The solution to the problem of oxide accumulation in the 
feed holes was more complex and required major changes in the 
polysilicon feed route through the growth system lids and shield 
stack. The new design is fully compatible with feed stock pellets 
prepared in the JPL shot tower developed by Kayex. 

These improvements permit melt replenishment to be interrupted 
at will and growth runs to be carried over to the next day, or several 
days, without oxide buildup. This in turn increases the amount of 
growth data which can be obtained in a single run and produces a saving 
in crucible and silicon costs. 



4. SUMMARY, CONCLUSIONS, AND FUTURE WORK 


4.1 Summary and Conclusions 

The thermal stress models have been successfully applied to the 
design of new low-stress growth configurations. The J419 and J460 
configurations represent the first two steps in achieving substantial 
reductions in buckling stresses, a factor of three overall, and have 
yielded substantial increases in the width of high-quality web crystals 
which have been grown. Growth experiments can now be carried out under 
controlled and constant growth conditions with the ability to maintain 
the melt at the desired level through continuous replenishment. The 
combination of thermal modeling with carefully controlled experiments 
can be expected to produce further substantial increases in the area 
throughput of web crystal of the quality required for high-efficiency 
solar cells. 

4.2 Future Work 

Thermal modeling combined with experimental evaluation will 
continue to be applied to the design of lower stress, higher throughput 
growth configurations. In addition, we will be studying dynamic thermal 
trimming as a method for optimizing both the requirements for growth 
during the initial transient Involving the start of growth and widening 
to desired width, and for growth at steady-state for increased area 
throughput. Initial efforts will be directed toward changes in melt 
level and the position of the top shields during web growth. Both 
parameters affect both stress and growth velocity, as well as the 
quality of crystal start. Thus, the objective would be to change the 
configuration during growth to match the thermal requirements at each 
stage. Additional configurational parameters will be examined as 
results and experience dictate. 
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5. NEW TECHNOLOGY 


No new technology is reportable during this period 
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APPENDIX I 

Web Temperature Computer Model 


1. INTRODUCTION 

Thermal stress and buckling modeling of dendritic web crystals 
requires temperature distribution data along the growing web. Because 
of the difficulty of measuring web temperature, we instead compute it 
from the furnace lid and shield geometry and temperature distribution. 
This computation is performed by a computer code called ’'RIBBON" which 
integrates the required heat transfer equation for the given radiative 
web environment. A description of this program and data input follows. 

2. CONSTRUCTION OF MODEL 


2 .1 Differential Equation 

The purpose of the RIBBON program is to determine the 
temperature, T, of the growing silicon web as a function of the 
distance, x, from the lower edge of the furnace lid. RIBBON 
accomplishes this task by integrating the heat conduction equation: 


P GpV 


dT 

dx 


d_ 

dx 


C a dT 'i 2q 
T dx b 


[A-l] 


where p * density 

Cp = specific heat 
V = web pull velocity 


a * 318 W/cm 
b = web thickness 

q = heat flux from one side of the web. 
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One of the major tasks of the RIBBON program is to calculate the 
geometric form factors for the term "q". A detailed discussion of this 
computation is found in the last quarterly report (DOE/J PL-955843/82/6) . 
For the purpose of recognizing the degree of nonlinearity of equation 
A-l, we note that q can be expressed in the form 

q * oET 4 + f ( x) , 

where a is the Stephan-Boltzman, E is the emissivity of silicon web, and 
f(x) is a function of position x on the web and the geometry and 
temperatures of the lid and shields of the furnace. Above 20 cm from 
the lid, the web enters a chimney and thus no longer "sees" the lid and 
shields. For a simple approximation, we say that it "sees" only the 
ambient temperature, T fl . Thus, for x > 20 cm, 

q - oe (T 4 - T^ 4 ) [A-2] 

The longer the ribbon grows, the closer its temperature approaches the 
ambient. This fact together with the initial condition that the web 
starts to grow at the silicon melt temperature, T m = 1685°K, gives us 
the boundary conditions 


T(x) = T 
' o m 

T(°°) = T a . [A-3] 

The position x * x Q is the growth front and has a value equal to the 
negative of the input parameter LIN in the program. 

Differential equations with boundary conditions are more 
difficult to solve than those with initial conditions; they generally 
must be solved iteratively. Since most numerical integration routines 
are written for first order equations, we transform the second order 
equation (A-l) into two first order equations with the substitutions: 
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u 2 - p Cp V T m /a 

H " eo ( ab ) 

TT(1) - T/T m 

TT(2) « || (1/T) [A-4] 

Equation A-l can be expressed now as a system of first order equations: 

TTP(l) - TT(i) • TT(2) 

TTP(2) =* u 2 TTP(l) + 3 2 Q 


where 

TTP(i) - — [TT(i) ] i = 1,2 

and 

Q = 2q/eo . [A-5] 

The nonlinearity of these equations makes them unstable as small 
errors in the steps of the integration are quickly magnified in the 
succeeding steps. We tried several different methods of numerical 
integration but found the simple fourth order Runge-Kutta method to be 
the most stable. Even so it is necessary to use double precision to 
obtain reasonable results. In practice, an initial slope of the 
temperature is guessed. If the slope leads to a curve which gives below 
the ambient temperature, then the slope is increased for the next 
guess. If the resulting curve goes above the silicon melt temperature, 
then the initial slope is decreased. This iteration continues until the 
double precision accuracy of the choice of initial slope is exhausted; 
in other words, there is no way to choose an initial slope between the 
high slope and the low slope since they are identical to 16 places. 
Generally, the integration curve does not blow up (or down) until 
20 - 30 cm. In this case the integration from 0 to 10 cm (the length of 
the buckling finite element model) is fairly accurate. If it blows up 
before this length, this may no longer be true. (A smaller integration 
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step size, HO, may improve this problem.) Also, if the integration 
remains between the melt and ambient temperature for a longer length, 
the integration may not be accurate because other values of the initial 
slope may lead to different temperature curves which also remain 
bounded. The boundary condition, T (») - i a> must then be applied. 

Thus, it becomes necessary to examine the asymptotic expansion of 
equations A-l and A-2 . While the numerical integration of equation A-l 
or equation A-4 starting with initial values at the melt interface will 
be accurate for small values of x, the asymptotic expansion can be 
expected to be accurate for large values of x. Hopefully, their regions 
of accuracy will overlap; if this is not the case for some problems, 
then some approximation to an intermediate solution might be required. 

2.2 Asymptotic Expansion 

For large x, it is most convenient to change the differential 
equations A-l and A-2 into a system of first order equations with the 
following substitutions: 


a * p CpV T a /a 
y - 2eo T„ V aba 2 

a 

y L =* T/T a 

y 2 * U/“« £ 


[A-6J 


Thus, 


dy^ 

dx 

^2 

dx 


a y x y 2 

. , 4 

a y i y 2 + TOt (y i 


1) 


[A - 7] 


If one attempts a formal power series solution of equation A-7 in 
negative powers of x (so that they are finite as x approaches infinity), 
one would obtain the trivial solution Y^ * 1 and Y 2 ■ 0. This solution. 
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however, does not give the general asymptotic expansion of equation 
A-7 . Theory (Wolfgang Wason, Asymptotic Expansions for Ordinary 
Differential Equations, Interscience, N.Y., 1965) shows that the general 
solution is a function of two parameters: C exp ax(l-/l+l6y)/2 and 

D exp ax(l+/l+16y)/2, where C and D are arbitrary constants. Since the 
second parameter becomes infinite as x approaches infinity, we let 
D vanish and look for functions of the first parameter. The substitution 

5 = exp [ax( 1-V 1+1 6y ) / 2 ] 


transforms the equation A~6 into 


y x ' 5 (1 - A+16y) - 2 Xi y 2 

y 2 'S (1 “ /I+ISy) = 2 y x y 2 + 2y (y^ - 1) [A-8] 


where the primes represent differentiation with respect to the variable 
We can now seek a power series solution of equation A-8 in the form 


y l = a 10 + a ll + a 12 + *•* [ A “ 9 ] 

2 

y 2 — a 2 i (C?) + a 22 **" ••• 


where a^Q = 1 from the boundary condition and a^ may also be chosen 
equal to unity since C is an arbitrary constant. Equating the 
coefficients of (CC) n on both sides of equation A-8, we obtain 


n a (1 - ✓l+l'6y) 
In 


n a„ (1 - A+16y) 
2n 


n-1 

2 y a a 

L n lm 2(n-m) 
m=0 


[A— 10 ] 


n-1 

2 l a, a 
m^O 


(n/2) -1 


+ 2y [b^/ +2 l b_ b ] 

m=0 


lm ~2(n-m) l "n/2 m n-m 
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for n even or 


na 2n 


n-l (n-l)/2 

2 L a lm a 2(n-m) + ^ J n b m b n-m 
m*0 m*u 


[A-10] 


for n odd, 


where 


b i ’ a l(i/2) 2 + 2 ^ a iJ a l(l-J) ’ 1 eVan 
(i-l)/2 

* .l Q a lj a Ki-j) , i odd 


b “ 1 
o 

From these equations, the a^ n and a 2 n are determined from the values of 
a lm and a 2 m , where m < n. In this way, we find the coefficients of the 
asymptotic expansion A-9. The radii of convergence of these power series 
may be found from the formulas: 

R = Lim | a n /a n+ i I if it exists 

and Lim sup I a n I = a = 1/R 

n + “ 

The next steps in the procedure are to choose a value of C and then a 
value of S > 20 cm. such that equation A-9 converges. From equations A-4, 
A-6, and A-9, initial conditions for equation A-5 may be found at the 
value of x corresponding to that of Equation A-5 can now be integrated 
backward to x * x Q using the same Runge-Kutta method as before. These 
steps are interated by increasing or decreasing the choice of C 
accordingly as T(x Q ) is less than or greater than T m . 
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3. COMPUTER PROGRAM INPUT/OUTPUT OPTIONS 


The input parameters of the RIBBON program are divided into two 
sets — the first to define the geometry and the second everything else. 

In this way, Calcomp plots may be obtained for any number of geometric 
configurations in each run. In each plot, up to three different graphs 
(black, green, and red) of the second derivative of aT can be obtained for 
different nongeometric parameters of web growth run.* A sample Calcomp 
plot is illustrated in Figure 2. Besides the aT second derivative curves, 
it diagrams the lid and shield geometry. 

The following input records are read in the first data set. The 
data are read in free field format with either spaces or commas separating 
them. Variables beginning with I-N are integers (except for LIN) and all 
others are double precision. 

RECORD 1: 

IGRAPH » Input 0 if there are to be no Calcomp graphs in the run, 
otherwise input 1. 

RECORD 2: 

NS * Number of geometric elements. These include the lid, the 
shields, and their separating gaps. 

JQC » Parameter to correspond to the type of WECAN element 

= 2 for quadratic elements 

= 3 for cubic elements 

JN * Number of WECAN elements in the x direction (generally has been 
set to 23). 

NE « Number of WECAN elements in the y direction (generally has been 5 
for two-dimensional elements and 7 for three-dimensional ones). 


*Here a stands for the thermal expansion of silicon and is represented 
by ALPO + ALPI * TEMP in the program; it is unrelated to the a in 
equations A6 and A7 . 
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NG * One plus the number of data sets for a given geometry. 

LIN * Distance (cm) of the growth front below the lower edge of the 

lid. 

EMAG 3 Magnification of the JNth element in the x direction relative to 
the first WECAN element at the growth front (generally has been 
8 ). 

ELL = Length of ribbon modeled on WECAN in centimeters (generally has 
been 10). 

RECORD 3: 

H(I),I = 1, NS = Height (cm) of ith horizontal surface above the lower 
surface of the lid. 

RECORD 4: 

Y(I), 1=0, ..., NS+1 = Half width (cm) of the central gap of the ith 

horizontal surface. 

The H and Y parameters are illustrated in Figure 3. This completes the 
first data set. 

The second data set consists of the following records: 

RECORD 1: 

IS =1 for forward integration 

= -1 for backward integration from asymptotic expansion (not yet 
implemented) . 

HO = integration step size (generally set to 0.01). 

ROUT = 6 for printed temperature data in WECAN stress input formalism 

= 7 for punched WECAN data. 

KQU = 0 for no plot of the otT second derivative. 

= 1 for a plot 

EPS = emissivity of silicon web. 
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* 


A ■ thermal conductivity of silicon multiplied by its temperature 

(318 W/cm). 

B - web thickness (cm). 

V ■ web pull velocity (cm/min) 

POMIN - minimum initial value of TT(2). If not known, set to -1. 

POMAX * maximum Initial value of TT(2). If not known, set to -0.01. 

CMIN ■ minimum value of the parameter in equation (A-9). If not known, 

try 1. 

CMAX =• maximum value of the parameter C. If not known, try 20. 

RECORD 2: 

TS(I),I * 1, ... , NS+1 = Temperature (degrees Kelvin) of the ith 

geometric element. 

These RECORDS 1 and 2 may be repeated up to a total of three sets. The 
data sets 1 and 2 can be repeated an indefinite number of times. 
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APPENDIX II 
Computer Code Listing 


BRUN, /RNPT JSS05,09F40RKADY20*SCHRUBEN*5*100/1500 
SHDG RIBBON . ’ 


FTN N $R?wi V 08/25/82 
1. C 


-i: 

6 . 

7. 

8 a. 

9. 

10 . 

11 . 

12 . 

13. 

14. 

15. 

]*: 


-20 


:50 

SILICON WEB - ANALYSIS OF HEAT LOSS OF MOVING WEB 
IMPLICIT DOUBLE PRECISION (A-H,0-Z> # 

DOUBLE PRECISION LIN *LAMB »TS <0 :2Q> ,H(0:20) * V(0 :20> , XS CO slOO) , 
1 X (0 :60) 

.DIMENSION TTC2) ,TS4 <20 > ,T< 100) , TP (100) , TPPC 100) , ATPPC100) ,Y2(t 
1 *TTP(2) 

INTEGER 
. INTEGER 
REAL 

THE ABOVE 
COMMON /CNI/ 

.COMMON ' 

DATA 


IXS(40).ISR(0:60.20> 

ER TITLE (1 1 ) s LAB ELX ( 2) ,LABELY(2 ) * I COLOR (4 ) 
XRC150) .VRC4.150) .XLEN 
BOVE ARE GRAPH PARAMETERS 

U2*BETA2.ISR,TS4,V2,H,ICHIN,NS,NS1 


18. 

19. 

20 . 
21 . 

if: 

24.. 

25. 


28. 

it: 

31. 

■H: 

34. 

35. 

36. 

37. 

38. 

39. 

40. 

41. 

»: 

44. _ 

45. 

«: 
48. 
4 9. 




52. 
5 3. 

54. 

55. 


II: 


58. 

59. 

60. 
61. 
62. 

63. 

64. 

65. 

66 . 
6 7. 
68 . 


N /CNJ / TA .1 FLAG , „ 
TITLE/'SECO*,'ND D','ERIV','ATIV V 
'EMPE * RATU * RE*5 , 00*/ 
LABELX/'X-AX*,'IS */ 

LA9ELY/ V— AX . IS. / .. 


’E OF',' ALP','HA*T', 


C, 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c. 

c 

c 

c 

c. 

c 

c 

c 

c 

c 

c 


-L 

WE 


DATA 

DATA I COLOR /I .1.2*3/ 

KI* 5 
KO* 6 

IGRAPH * 0 FOR NO CALCOMP ..... . 

* 1 FOR SOME CALCOMP 

NG IS ONE PLUS NUMBER OF DATA SETS FOR A GIVEN GEOMETRY (4 MA) 
NS IS NUMBER OF SHIELDS PLUS LID SECTIONS PLUS GAPS INBETWEEN 

. TSClS IS THE TEMPERATURE OF THE ITH SHIELD OR GAP 

HCI) IS THE HEIGHT OF THE ITH HORIZONTAL SURFACE ABOVE THE LOl 
SURFACE OF THE LID 

Y l I J IS THE DISTANCE OF SHIELD EDGE TO WEB 
. (YCO) -IS HALF-WIDTH OF MELT SURFACE* Y(1> IS HALF GAP OF LID*. 
YCNS) IS DISTANCE OF INNER EDGE OF TOP SHIELD TO WEB. 

V (NS ♦ 1 I IS DISTANCE OF OUTER EDGE OF TOP SHIELD TO WEB) 

JQC *2 FOR QUADRATIC ELEMENTS 

=3 FOR CUBIC ELEMENTS .. . . 

IS = 1 FOR FORWARD INTEGRATION 
*-1 FOR BACKWARD INTEGRATION 
JN IS NUMBER OF ELEMENTS IN A ROW ALONG THE WEB AXIS 
jLIN^I S DJgJANCE 0F GROWTH FRONT. BELOW LID . 

EMAG IS MAGNIFICATION OF LAST ELEMENT RELATIVE TO FIRST NEXT 1 
ELL IS LENGTH OF WEB TO BE MODELED 

ELL <20 _ 

FOR 2D ELEMENTS 
FOR 3D ELEMENTS 
1685 .DO 

1685. DO ..... _ 

1./TM 


ASSUME 
NE* 5 
NE* 
TSCO 
TN* 
TMR 


U 


LAMB* 1804. DO 
TM4* TM**4 
CP* 0 .981 1 DO _ 
RHO* 2.30D0 ~ 

SIGMA* 5. 670— 12 
HCO)* 0 
X<0>* -1 
ISR(O.O)* 0 
XO* X(Q) 

ISR(0,1>= 1 
XS(O)* X0„ 

RE AD (K I **) IGRAPH 
IF (IGRAPH .GT. 0) 
VR (1 .1 )* -1 
XR ( 1 5 * -1 
XR (2i* -1 

DO 998 ISUPER*1 *10 
READ (KI ,*,END=999) 


CALL PLOTINd .*1 •) 


READ(KI .*,END=999> NS, JQC, JN,NE,NG,LIN,EMAG,ELL 
WRITE («. S) , JQC * JN.NE *NG. LIN * EMAG * ELL 

FORMAT (*U NS**, 14,* JQC* ,14,* JN* ,14, NE*' t I4,' NG*', 

1 . ' LIN*',E12.5,' ENA 6*', El 2.5* *LL*%E12.5> 


14, 


NS1* NS* 1 
H ( N S 1 ) * 20 
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V» v/*v/»uiviv»u**'.ra>.04CNwrow *“04 rvirvjr\jfor\j»\j-» 


RIBBON 


69. 

70. 

71. 
_7H. 

73. 

74. 

75. 

,. 76. 

77. 

78. 

79. 

80. 
81. 
82. 

83. 

84. 
65. 
86 . 

87. 

88 . 

89. 

90. 

91. 

. 92. „ 

93. 

94. 

95. 

96. 

97. 

98. 

99. 
1Q0. 
101 . 

181 : 
104. ... 
IQS. 

W: 
108. „ 
109. 

HI: 

112 . 

113. 

114. 

115. 

1 1 6. 

117. 

118. 
119. 

1 20.„ 

121 . 

122 . 

123. 

124. 

"T 2 5 • 

126. 

127. 

128. 

129. 

130. 

131. 

132. 

133. 

134. 

135. 

136. 

137. 

ill: 

140. 

141. 

Hi: 

144. 


7 

8 

10 


20 


READ(KI.*) CH(I).Ix1 f NS) 

SHIHijti/I JlHiiSjISH, 

FORHAT<*CY» \/,<D14.5>> 

00 CO I* 0.NS1 
Y2 < I J * Y<I)**2 
CONTINUE 

YR(1,2>* SNGL(YCO)) 

ALPO* 2.8457E-6 
ALP1* 9.796E-1Q 
ALP2* 2*ALP1 
X * 0 

00 SO 1= 2, NS 
II* 1-1 
XLS* XO 

XLS IS THE LOWER BOUND OF TS(I> RADIATION ON WEB 
HI* H % I » 

YI* V ( 1) 

DO 20 4* 1.11 

IF < Y I .GT. Y ( J ) ) THEN 

TEMP* (H(4)*YI-Y (4)*HI)/(VI-Y<4)> 

. IF (TEMP .GT. XLS ) XLS= TEMP 

END IF 
TONTINUE 

IF (XLS .LE. XO) THEN 
ISRCO.I)* 1 


ELSE 


ISR(O.I)* -1 
DO 25 4=K*0»-1 
41* 4+1 

IF (XLS .GE. XS(J>) 
XS(J1)= XLS 
IXSCJ1>= I 
GO TQ 26 


THEN 


ELv'C 


25 

26 


END 

XUS 


IF 

ZQ 


END IF 
CONTINUE 
K* K+1 


XS<41)= XSCJ) 
IXSCJ1)* IXS(J) 


XUS IS THE 
i-2 ) 
VI* Y ( I 1 ) 


HI* HCI- 


UPPER BOUND OF TS(I1) RADIATION ON WEB 


30 


DO 30 J= I, NS . 

J1* J-1 

IF ( Y I .GT. V ( J ) ) THEN 

TEMP* CYI*H(J1)-HI*Y(J))/(YI-YCJ) ) 
If.. (TEMP .LT. XUS) XUS* TEMP 

END IF 
CONTINUE 

IF (XUS .LT. 20) THEN 
DO J5 J= K,0,-1 


J1* 


(XUS ,6E. X$(J)) 
XS 1 41 ) * XUS 
IXSC41)* II 
GO TO 36 


THEN 


35 

36 

50 


60 


ELSE 

XS(41)= XS (4) 

IXS (41 )* IXS (4 ) 

END IF 
CONTINUE 
X* K + 1 

END IF 
CONTINUE 
IB* 0 
I* 0 

DO 70 L* 1 » K 

IF (XS(L)-X(I) .LT. 1.D-4) THEN 
IXSL* I X S ( L ) 

ISR(I.IXSL)* -ISR(I.IXSL) 
. SO T0..7Q . 
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RIBBON 


145. 

146. 

147. 

148. 

149. 

150. 

151. 

152. 

153. 

154. 

’?!: . 

157. 

ill: 
160. . 
161. 
162. 

163. 

164. 

165. 

166. 

167. 

168. .. 

169. 

170. 

171. 

172. 
17,3. 

174. 

175. 

176. 

177. 

178. 

179. 

180. 
181. 
182. 

183. 

184. 

185. 

186. 

187. 

188. 

189. 

190. 

191. 

192. _ 

193. 

194. 

195. 

196. 

197. 

198. 

199. 

200 . . 
201 . 

m: 

204. 

205. 

206. 
207. 
208< 

209. 

210 . 
211 . 
212 . 

213. 

214. 

215. 

Ii$: 

218. 

219. 

220 . 
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END IF 
II* I 
I* 1*1 

DO 65 J* 0 , NS 

ISR ( I » J ) = ISR(I1tJ> 

IF CH(IB)-XSCL) »GT • 1.D-4) THEN 

x c i > * xs cl) 

IXSL= IXS(L> 

ISRCI , I X SL ) * -ISRCI, IXSL) 


ELSE 


END IF 
70 CONTINUE 
DO 80 L« 

IF (HCL) 
II* I 

.. 1= 1*1 

X(I)= HCL) 


XCI )* H C IB) 
IB* IB*1 ... . 

ISRCI *0 ) * IB 
GO TO 60 


IB t NS 

.LT. X Cl) ) 


GO TO 81 


. DO 75 J* 1.NS 
75 ISRCI, J)* I SR (II, J) 


THEN 


__ 80 ISRCI, 0)* L* 1 .. 

81 CONTINUE 
LAS* I 

IF (XCI) .NE. 20) 

LAS* 1*1 . .. . 

X (LAS) = 20 
END IF 
LAS* LAS* 1 

X (LAS ) * 30 

WRITE (KO, 1000 ) (XCJ) ,J*0,LAS> 

WRITE (KO ,1001 ) CXSCJ 5 , 1 XS (J > .J*1 ,K) 

WRITE CK0,1 002 ) CClSR (L, J > , J *0 ,NS ) , 1*0 , I > 

DO 90 I* 1 • LAS 
K= I 

IF((XCI)*LIN) .GT . 1.D-8) GO TO 91 

90 CONTINUE 

91 CONTINUE 

KO* K 

RA= ENAG**(».DO /CJN-1)) 

EL* ELL*C i-RA)/C1-RA**JN) 

JNN* JQC*JN+1 . 

KK* 2 

DO 110 KK I * 0 ,K 
KKKS* KK 1*1 
KKK* KKI 

I F ( (H ( KK I )♦ LIN )• GE . 0. ) 60 TO 111 

KK* KK ♦I 

x, :kk)= snglchckki)) 

_ YHC1,KK>* SNGLCYCKKD) ... ... 

KK* KK*1 

XRCKK)* SNGLCHCKKI)) 

YR (1 ,KK)= SNGLCY(KKKS) ) 

110 CONTINUE 

111 KKS* KK 
KKS1 * KK*1 

DO 997 IC* 2 ,NG 

READ C K I , * ) IS,HO,KOUT,KOU,EPS,A,B,V,POHIN,POMAX.CMIN,CM 

WRITE CKO, 82) IS ,H0 ,KOUT .KOU, EPS , A ,B ,V ,POMIN ,POMAX ,CHIN , 
C KOUT = 6 FOR PRINTED STRESS DATA 

C * 7 FOR PUNCHED STRESS DATA 

C KOU > Q FOR GRAPH 

82 FORMAT C'O IS*', 13." H0*',E12.5, 

1 ' KOUT*', 13," KOU*', 13,' EPS* ',E12.5.' A* ',E1 

2 E12*5, V* ',E1 2 .5 , / • ' POMIN* '.D25.18,' POHAX* 

_ .. .3 # /,* CMIN* ,[>2 5 .18 ,' CHAX* , D25.18) 

READCKI,*) (TSCI) .1=1, NS1 ) 

WRITE ( KO ,6) CTSCI) ,I*1,NS1) 

6 FORMATC'OTS*',/, (014.4)) 

DO 83 I* 0.NS1 
TSCI)* TS(I)*TMR 
TS4CI)* TS C I ) **4 

83 CONTINUE 

IS1* M0DC1-IS.4CC) 
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RIBBON 


60 TO 102 


► » + l 

iW: 


NO* HO*IS 
TA* TS(NSI) 

V* V/60 

BE.TA2*. EPS*S16NA»TH4/(A*D> .. „ 

BETA* DSQRT(2*BETA2) 

U2* RHO*CP*V*TN/A 
DO 101 J* 1,500 

I FLAG* 0 ... . _ . . . 

K* KO 
OX* EL/ JOC 

IF (IS .EG. -1) 60 TO 400 

TT (1 ) * 1 .. ... 

XX* -LIN 

POOO* >500* (POHAX+PONIN) 

IF (POOO .EG. POO .AND. 4 *NE. 1) 

— ???I>S°°8o - 

KK* KKS + 1 
KKK* KKKS-1 

11*1. 

12* JNN-1 
X2* -LIN 
ICHIN* 1 

60 TO 409 

400 CONTINUE 

409 CONTINUE 
TENP* TT(1)*TN 

..... TENPP* TT (2 ) *T ENP ... . __ 

DO 450 1*11 ,12, IS 
XI* X2 
X2* X1+DX 

IF (I .EG. 2 .AND. IS .EG* rl > X2* X(1) 

KFLAG* 0 
T(I)» TENP 
TP(I)* TENPP 
JISCI)* XI 

IF (IS*(X(K>-X2)) 410,411,412 

410 CALL FCNJ(X1.X(K) ,TT,TTP,HQ,K) 

IF (IFLA6 .NE. 0> 60 TO 600 

XR (KK) * SN6LCX1) 

TENPPP* TENP*TTP(2)*TENPP*TENPP/TENP 
ATfcfel* (ALPQ+ALP2*TENP)*TENPPP*ALP2*TENPP*TENPP 
TENP* TT(1)*TN 
TENPP* TT (2>*TENP. 

YR(IC,KK>* ANAXM-1..ANIN1 (8 . ,SN6L (ATEN ) *500 • ) ) 
IF (KFLAG .EG. 0) THEN 
TPP ( I ) = TENPPP 

.. ATPP(I)* ATEN 

END IF 

IF (H (KKK) .EG. XI) THEN 
KKK* KKK + 1 

_ KK* KK+1 

XR (KK ) * SN6L ( XI ) 

TR(IC,KK)* YR (IC ,KK— 1 ) 

END IF 

KK* KK+1 

” KFLAG* 1 ' 

XI* X (K) 

411 K* K+I S 

IF ( IS* (X(K)-X2) ) 410,411,412 

412 CALL FCN4(X1,X2.TT,TTP,H0.K) 

IF (IFLA6 .NE. 0) 60 TO 600 
TENPPP* TENP*TTP(2)+TENPP*TENPP/TENP 

_ ATEN* (ALP0+ALP2*TENP)*TENPPP*ALP2*TENPP*TENPP 
" TENP* TT ( 1 ) *TN 
TENPP* TT(2)*TENP 
XR (KK ) * SN6L (X 1 ) 

YR(IC,KK)* ANAXK-1..ANIN1 (8. ,SNGL( ATEN) *500. ) ) 


IF 


END 

KK* 

IF 


(H (KKK) .EG. Xl) THEN 
KKK* KKK+1 
KK* KKil 

XR (KK )= SN6L( XI > 
YR(IC.KK)* YR (IC ,KK— 1 ) 
IF 

KK + 1 


(KFLAG 


0) 60 TO 450 
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RIBBON 


mt I 

i\ 


450 


TPP f I ) - TENPPP 
ATPP(I)= ATEN 

IF (MOD(I.JQC) .EQ • IS 1 > DX* DX*RA**IS 

.CALL FCNI(X2,TT,TTP,K) ...... _ . 

13* I2 + IS 
T (13) * TT(1)*TN 
I* 13 
XS(I)* X2 
TP (I ) * TT (2 ) *T (I ) 

TPP (I ) * TCI)*CTTC2)**24TTPC2)) 

ATpptn* (alpo+alp2*t\ i r 
.. XR(KK) = SNGLCX2) . . 

YR(IC,KK)= AMAX1C-1 • » AMINl (8 . ,SNGL (ATPP C I ) ) *500 . ) ) 
IF (IS .EQ. -1) GO TO 600 
“ 500 JX* 


>)*TPP(I)*ALP2*TP(I)**2 


... K t LAS 

IF (4K »EQ. LAS) ICHIH* 0 
CALL FCNJ(X2,X(JK),TT,TTP.H0,JK) 
IF (IFLAG .NE. 0) GO TO 600 
X2* X ( JK ) 

— 500 CONTINUE . .... - . . 

600 CONTINUE 

780 WRITE(K0,1002) IFLAG 

WRITE (K0,1005) XI ,X2,TT,P00,C 
— 100S F0RHATC3D2S.18) 

IF (IS .EQ. -1) GO TO 95 
IF (IFLAG) 93,92,94 
92 GO TO T02 
- 93 POHIN* POO 


94 


GC TO 101 
POMAX* POO 


849 

8490 

852 


GO TO 101 

-.95 CONTINUE ...... 

101 CONTINUE 

102 VV* -60*A* POO/(RHO*LAHB) 

WRITE <K0. 7100) VV 

7100 F0RMAT('1',5X,'VV= '.E14.7) 

DO 849 1=1,5 

YIELD* 2.57E-11*EXP(494597T(I>) 

WRITE (KQ, 8490) YIELD 

FORPATC* YIELD* ,E15. 8) ..._ . _ 

DO 852 1=1, JNN 

WRITE(KO,1000)XS(I),T(I),TP(I),TPPCI),ATPP(I) 

DO 850 1*1, JN 
II* JQC*(1-1)+1 
“ T ( 1 1 ) 

T ( II+1 ) 

T(II ♦?) 

T(II43) 

800 4*1, NE 
II* J+CI-1 ) *NE 
IF (JQC .EQ. 3) GO TO 795 
WRITE (KOUT, 1006) II.T1 ,T1 ,T3,T3,T1 
IF (NE .EQ. 5) THEN 

WRITE (KOUT, 1008) T2,T3,T2 

ELSE 

. _ WRITECK0UT.1007) T1,T3,T3,T1,T2,T3,T2,T1,T1,T3 

WRITE (KOUT, 10 08) T3 ,T1 ,T2 ,T3 ,T2 

END IF 
GO TO 800 

WRITE (KOUT, 1006) II,T1 ,T1,T4,T4,T1 
WRITE(K0UT,1007) T2,T4,T3,T1 ,T3 
WRITE (KOUT, 1008) T4,T2 
CONTINUE 
CONTINUE 
YRCIC .1)= 8 
DO 860 KKI* 

YRCIC, KKI)» 

CONTINUE _ 

CONTINUE 
ISW* -1 
K= KKKS-1 
DO 865 KKI* 

IF (XR(KKI) .NE. SNGL(H(K))) THEN 


795 


860 

997 


2.KKS 

YR (IC ,KKS + 1 ) 


KKS1 ,KK 

.NE. SNGLCH (K ) ) ) 
YR (1 ,KKI )* YR (1 , KK 1-1 ) 

K* Kt H AXu (0 t ISW ) 
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ELSE 
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RIBBON FCNI 


R AYS*R 1B80N (1 J.FCNI <0 > 


1 

2 

3 

4 

5 

6 

7 

8 
9 

10 

11 

12 

13 

14 
1 5 
16 

17 

18 

19 

20 
21 
22 

23 

24 

25 

26 

27 

28 

29 

30 

31 

32 


50 


60 

100 

1000 


SUBROUTINE F C NI (X , TT , TTP , K > 

IMPLICIT DOUBLE PRECISION (A-H.O-Z) 

DIMENSION TT(2),TTP(2) ,TS4 <20 ).Y2 (0:20) ,ISR (0:60,20) ,H(0 
COMMON /CN1/ U2,EETAZ,ISR,TS4,y2,H,ICHIH,NS,NS1 * ’ U 

TTP (1 ) = TT (2 ) *TT ( 1 ) 

Q * 2*<TT(1 )**4-TS4(NS1) > 

IF (ICHIM .EG. 0) 60 TO 100 
0= Q+TS4 (NS1 ) “1 
T4= 1 

N= I S R (K 1 ,0) 

N1= N* 1 

DO 50 KK =1,N 

IF (ISROC1.KK) .LT. 0) 60 TO 50 
XX = HCKK^Ii^X 

G= Q* <TS4(KK)-T4>*XX/DSQRT(XX*XX+Y2(KK)> 

T4= TS4 (KK ) 


: 2 


CONTINUE 
XX = H (N ) -X 
Y2K= Y2(N) 
DC 60 KK= 


N 1 , NS 


mw uu ms- niyno 

IF (ISR(KI.KK) .LT. 0) GO TO 60 

Q* Q+ C TS4 CKK > -T4 ) *XX/ DSQRT (XX* XX* Y2 K ) 

X X= H (KK) — X 

Y2K- Y2(KK) 

T4= TS4 (KK ) 

CONTINUE 

Q= Q ♦ (TS4(NS1)-T4)*XX/DSGRT(XX*XX*Y2K) 

TTP (2 ) = U2*TTP(1 )+BETA2*Q 

RETURN 

FORMAT (8E1 5.8) 

END 


SRJEFIN : 


R JESPF.RJEFIN . . . 


8SS-PPF-SS.RJESFIN 



RIBBON FCNJ 


RAYS*RXBBON(1).FCNJ(0) 


1 


4 

5 

6 

7 

8 
9 

10 

11 

12 

13 

14 

15 
1 6 

17 

18 

19 

20 
21 

if 

24 

25 

26 

27 

28 

29 

30 

31 

32 

33 

34 

35 

36 

37 

38 

39 

40 

41 

42 

43 

44 

45 

46 

47 

48 

49 

50 

51 

52 

53 

54 


SUBROUTINE FCNJCZI.ZZ.Y.YP.HO.K) 
IMPLICIT DOUBLE PRECISION (A-H.O-Z) 
DIMENSION Y (2) ,YT (2) ,PY(2 ) , YP(2 ) 
COMMON /CN J / TA , I FLAG 
2 - Z1 
H= HO 

NI s DMAX1 (1 . I DINT < (Z2-ZD/H) ) 

DO 500 1= 1 f NI 
IF (I .EQ. NI) H= Z2-Z 
HH= H* • 5D0 


CALL FCNI CZ.Y.PY.X) 
IF (1 .EQ. 1) THEN 


YPC1 > = PY (1 ) 

Y PC2 ) s PY ( 2 ) 

END IF 

RK 1 = H + PY (1 ) 

RL1= H+ PY (2 ) 

YT Cl )= Y(1 )+RKl+.5D0 
YT ( 2>= Y(2)+RL1*.5D0 
Z= Z+HH 

CALL FCNI (Z .YT.PY.K) 

RK 2= H*PY(1) 

RL2= H*PY(2) 

YT Cl )=Y <1 )+RK2*.5D0 
YT C 2 ) = YC2) +RL2+.SD0 
CALL FCNICZ* YT.PY.K) 

RK 3 * H + PY Cl ) 

RL3= H* PY ( 2 ) 

Y T ( 1 ) = Y C 1 ) +RK3 
YT C 2 ) = Y(2)+RL3 
Z = Z+HH 

CALL FCNI CZ, YT.PY.K) 

RK 4 = H+PYC1) 

RL4= H+PYC2) 

Y ( 1 ) = YCD + CRK1+RK4 + 2 + CRK2 + RK3) ) /6 
Y(2)= YC2)+(RL1+RL4+2*CRL2+RL3) ) /6 
IF ( Y (1 ) .GT. 1) GO TO 10 
IF CY ( 1 ) .LT. TA) GO TO 11 
500 CONTINUE 
780 CONTINUE 
RETURN 

10 CONTINUE 

IF C I F LAG .EQ. 0) 

1 UR I TE (6 » 1 000 ) Z » Y 
1 FLAG= 1 
RETURN 

11 CONTINUE 

IF ( I FLAG .EQ. 0) 

1URITE C6 .1000) Z.Y 
I FLAG= -1 

Q P T IIPM 

1001 FORMAT (1018 ) 

1000 FORMATC8E15.8) 

END 


3 HDG RIBBON FCNI .L»1 
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